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Abstract. We discuss how to compute certified enclosures for the eigenval- 
ues of benchmark linear magnetohydrodynamics (MHD) operators in the plane 
slab and cylindrical pinch configurations. For the plane slab, our method re- 
lics upon the formulation of an eigenvalue problem associated to the Schur 
complement, leading to highly accurate upper bounds for the eigenvalue. For 
the cylindrical configuration, a direct application of this formulation is possi- 
ble, however, it cannot be rigourously justified. Therefore in this case we rely 
on a specialized technique based on a method proposed by Zimmermann and 
Mertins. In turns this technique is also applicable for finding accurate comple- 
mentary bounds in the case of the plane slab. We establish convergence rates 
for both approaches. 



1. Introduction 
Let 17 C M^. The linearized ideal MHD equation 

pd^at) + Km = m - eo otm = m 

for displacement vector ^ : — > M.'^ and force operator 

K(, = grad[7P(div ^) -I- (grad P)-^] + -[Bx curl(curl(B x ^)) - (curl B) x curl(B x ^)] 

M 

arises in applications from plasma confinement in thermonuclear fusion. The con- 
stants fi and 7 here denote the magnetic permeability and heat ratio. The smooth 
function p is the density, P > is the pressure and B the divergence-free magnetic 
field of the given equilibrium, satisfying /xgradP = (curlB) x B. 

In the study of this equation a fundamental role is played by the eigenvalue 
problem associated to K. The appropriate Hilbert space setting ensures that K 
has a self-adjoint realization. A considerable amount of research has been devoted 
to the formulation of a rigorous operator theoretic framework for K and to the 
structure of the spectrum, [17] . Particular attention has been payed to the plane 
slab (plasma layer) and the cylindrical (plasma pinch) configurations [9l [TOl ETJ [1] 
where K is reduced to a block ordinary differential operator matrix. A systematic 
description (analytical or numerical) of properties of the eigensolutions turns out to 
be difficult even for these, the simplest configurations. This is due to the presence 
of regions of essential spectrum near the bottom end of the spectrum. 

The plane slab configuration has been the subject of thorough analytical inves- 
tigation and has become a benchmark model for the top dominant class of block 
operator matrix, see |26] and references therein. Precise eigenvalue asymptotics 
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can be found in this case by means of the WKB method 17,, §7.5] or by means 
of speciahsed variational principles, see [26l Theorem 3.1.4] and references therein. 
The cylindrical configuration is more involved due to the presence of singularities in 
the coefficients of the differential expression; however, eigenvalue asymptotics are 
known in this case, j20) . 

Specialized variational approaches are extremely useful for examining analytic 
asymptotics for the eigenvalues in the case of the plasma layer configuration. Un- 
fortunately, as they usually involve a triple variation formulation, it is arguable 
whether they are well suited for direct numerical implementations. 

If a sequence of subspaces is guaranteed not to produce spectral pollution, then 
the standard Galerkin method can be used. A prescribed recipe for avoiding spuri- 
ous modes when these subspaces are generated by the finite element method dates 
back to [221 H^- In this classical approach convergence is guaranteed, however, it 
is never clear whether a computed eigenvalue is on the left or on the right of the 
exact eigenvalue. In this respect the method is not certified. 

A technique for finding certified enclosures for the eigenvalues of K in the case of 
the plane slab configuration was considered in 24j based on the method proposed in 
[H [231 US] • The approach was based on computing the so-called second order spec- 
trum of K for given finite dimensional subspaces generated by the spectral basis. 
In the present paper we consider a further computational strategy which improves 
upon this technique in terms of accuracy. Our main approach is to combine two 
complementary Galerkin-type methods for computing eigenvalue enclosures which, 
by construction, never produce spectral pollution. 

For the plane slab, our method relies on the formulation of an eigenvalue problem 
associated to the Schur complement, this leads to highly accurate upper bounds. 
For the cylindrical configuration, a direct application of this formulation is possible, 
however, it cannot be rigorously justified. Therefore in this case we rely on a 
specialized technique based on a method proposed by Zimmermann and Mertins 
[l9] as described by Davies and Plum in [6l Section 6] . This approach is intimately 
related to classical methods, see [TTl [THl [Hj- We also apply this technique to the 
Schur complement and find accurate complementary lower bounds for the plane 
slab. 

In Section 2 we give a mathematical formulation of the MHD operators under 
investigation, and some of their spectral properties. In Section 3 we examine the 
approximation technique due to Zimmermann and Mertins. We present a formu- 
lation of this technique in terms of the Galerkin method which establishes both 
approximation and, importantly, the convergence of the method. Our main results 
are contained in Section 4. We present a highly efficient method for obtaining upper 
bounds for eigenvalues above the essential spectrum of top-dominant block opera- 
tor matrices, an example of which is the matrix K associated to the plasma layer 
configuration. We show in Theorem 14 . 41 that the convergence rate for this approach 
is the same as that achieved by the Galerkin method when applicable (below the 
essential spectrum). Our method is therefore extremely efficient. We also combine 
this approach with the Zimmermann and Mertins technique to obtain complemen- 
tary lower bounds for the eigenvalues. In Theorem 14.71 we use our results from 
Section 3 to obtain convergence rates for these lower bounds. In Sections 5 and 6 
we apply our results to the plasma and cylindrical configurations, respectively. 
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2. One-dimensional MHD operators 

The reduction process for the force operator, the precise constraints on the 
equilibrium quantities and the boundary conditions on fl, which yield the one- 
dimensional boundary value problems associated to the plane slab and cylindrical 
configurations, are described in detail in 17, §7.2 and §8.2], respectively. Through 
this reduction K becomes similar to a superposition of operators which are self- 
adjoint extensions of block matrix differential operators of the form 



(1) Mo - 



A B 
B* D 



acting on i^-spaces of a one-dimensional component. 

For the plasma layer, the components of A/q are explicitly given by 

' A = -p^'d^Pn{vl+vl)d^ + k^vl 

B = {{~ipa^dxPo{vl + vl) + ig)k^, {-ipg^d^pav^^ + ig)kii) , 

k'^vl + k'^v'^ /c^fcj|W^^ 



(2) 



D = 



kikiiv"^ kf.vl 



where 



Dom(A) = H^{{Q, l);podx) n HI{{0, I); p^Ax) 
Dom(B) - [i/i((0,l);poda;)]2 
Dom(i?) = [L2((0,l);poda;)]2. 

We assume that the Alfven speed Va, the sound speed Vs, and the coordinates of 
the wave vector k± and fcy , are bounded differentiable functions. We also assume 
that po and Vs are bounded away from 0. Following standard notation in the 
literature k'^{x) = k'^{x) + fcjj(a;) and g is the gravitation constant. In this case [26l 
Proposition 3.1.2] 

(3) Mo : [Dom(^) nDom(B*)] x Dom(S) — > [L'^{{0,1)] podx)f 

is essentially self-adjoint. Denote by Alg the closure of Mq. Then Ms is bounded 
from below, and the essential spectrum is given by the range of the Alfven frequency 
DqA:|| and the mean frequency v'^v'^k\\/{v'^+v'^), see jTTl §7.6]. The discrete spectrum 
always accumulates at -l-oo. We show in Example 1 2 . 1 1 that endpoints of the essential 
spectrum can also be points of accumulation. 

Example 2.1. Let k\\ = p = Va=Vs— g~l and k± = 0. Then SpeC(,ss(Ms) = 
{1/2, 1} and SpeCdi,(Ms) = {X^j^^, where 



\+ = l + k^n^ + VT+¥^r^+k^ =0{k^) and 

Ar = 1 + fc^TT^ - Vl + fc^TT^ + fc^TT^ - as fe ^ CX). 

Both A^ are positive and increasing in k. Also 



Range / dE\ ~ span{0j, }, Range / dEx = span{(/)^} 

J(-oo,3/4] J(3/4,3/2] 

and Range / dEx = span{(/)jj'}, 

J(3/2,oo) 
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where 



(4) 



sin kirx \ 1 f ^ 

— (kn cos knx + sin knx) I 



--|- (fcTT cos knx + sin kirx) 



a-A- -/ V 

(sin kiTX 

the constants chosen so that are normahzed. Note that {0^,0^,0^} 

is an orthonormal basis of [i^(0, 1)]"^. 

Example 2.2. Let po = I, kj^ = 1 , k\\ = 1, g = 1, Wa(a;) = ^7/8 — a;/2, and 
fs (a;) = -^/l/S + a;/2. The essential spectrum of Mg is given by 

Ran(t;2fc||) U Ran{vy,k^/{vl + v^)) = [7/64, 1/4] U [3/8, 7/8]. 

Below we will use the fact that d = niaxSpcc(Z?) = 1 + ■\/l7/32 < tt^. 

The plasma pinch configuration yields a differential operator Mq with singular 
coefficients. 



B = {{-idriP + lP)m^ + 2ib^k'-^, -idrlPk^) 

D = 



(5) 

m4,k^-fP kl^fP ) 

acting on [L'^ {{Q , Rq); r Ar)]^ . Here 9* = ^^.r, 

B = (0, b{r) sin{(j){r)), b{r) cos (f>{r)), 

P'{r) = -b{r)h\r) = -b{rf sin^ 0(r), 
r 

h{r) and 0(r) are smooth functions with 6'(0) = (/)(0) = 0, 

m . TO . 

kch — k cos (i) H sin and m,* = — cos (j> k sin 0. 

r r 

The indices i?o^ and m are integer numbers corresponding to the Fourier mode 
decomposition of K. 

In order to define rigourously the domain of Mq for this configuration, a fur- 
ther change of variables r = e'' is usually implemented, [51. Under this change 
of variables, Mq becomes similar to an operator acting on [L'^{{—oo,\ogRQ);dx)]'^ 
which is essentially self-adjoint in the space of rapidly decreasing functions at — oo 
vanishing at logi?o, [9^ Theorem 2.3]. Operator Mq is essentially self-adjoint in the 
pre-image of this space under the similarity transformation. We denote by M^. the 
unique self-adjoint extension of Mq in the latter domain. 

The original formulation ([S]) is numerically more stable for the treatment of the 
eigenvalues via a projection method. The finite element space generated by Hermite 
elements of order 3, 4 and 5, subject to Dirichlet boundary conditions at and Rq, 
considered below are C-'^-conforming and hence are all contained in Dom(Afc) for 
the benchmark equilibrium quantities considered in our examples. 
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The essential spectrum of Mc consists of an Alfven band determined by the range 
of fe^fc^, and a slow magnetosonic band determined by the range of {iP'k'^'-^P) / [li^ + 
7-P)j El Theorem 3.5]. As in the previous configuration, these bands are located 
near the bottom of the spectrum and +00 is always an accumulation point if the 
discrete spectrum. 

Example 2.3. Let F = 0, = 1, (/) ee 0, and fc = to = 1. Then SpeCc,g(Af) = {0, 1} 
(the point is the slow magnetosonic spectrum and the point 1 is the Alfven 
spectrum). On the other hand SpeC(jig(M) = {E"^ + 1 : J'i{E) = 0} where Jv{x) is 
the Bessel function of index v. 

3. Pollution-free bounds for eigenvalues 

The essential spectrum of both operators Mg and Mc is non-negative. Therefore 
unstable spectrum can only occur in the discrete spectrum. The eigenvalues below 
the bottom of the essential spectrum can be computed using the standard Galerkin 
method. Hence, the stability of the configuration can be determined by means of 
the Rayleigh-Ritz variational principle. 

By contrast, computing the eigenvalues above the essential spectrum is prob- 
lematic due to the possibility of variational collapse. The technique described in 
this section avoids spectral pollution and can be implemented on the finite element 
method. It gives certified enclosures up to machine precision for eigenvalues above 
the essential spectrum. In Section|3]we argue that this technique should be applied, 
not only to Ms, but also to its Schur complement. In order to keep a neat notation, 
we formulate the general procedure for a generic semi-bounded self-adjoint operator 
T acting on a Hilbert space %. 

3.1. Basic notation. Let the dense subspace Dom(r) C H be the domain of 
T. For ^ — minSpec(r), let t be the close bilinear form induced by the non- 
negative operator T — ji with domain Dom(t) = Dom(|r| 2 ). We denote the inner 
products and norms that render Dom(T) and Dom(t) with a Hilbert space structure, 
respectively by (m, v)t = {Tu, Tv) + {u, v), {u, v)t = t{u, v) + {u, v), \\ ■ \\t and || • 
Let £ he a subspace of Dom(T). For another subspace C, we denote 

s{£, c) = sup dist[(/), c] a Can 

0e£, 11011=1 

Si{£,£)^ sup disU[(j),C] if£cDom(i) 

4>e£, ||0||t=i 

St{S, C) = sup distT[0, C] if £ C Doni(r). 

0G£, ||0||t = 1 

Here and elsewhere dist, [(/),£] refers to the Haussdorff distance in the norm || • ||, 
between {(f)} and C. 

Below we establish spectral approximation results by following the classical 
framework of [5]. These results will be formulated in a general context for se- 
quences of subspaces £„ C H which are dense as n — >■ 00 in the following precise 
senses. We will say 

(£„) e A = A(/) dist[ii, £„] ^0 Muen 

(£„) e A(t) distt[M, £„] ^0 Vu G Dom(t) 

(£„) e A(T) ^ d\stT[u,Cn] ^0 Vw e Dom(T). 
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Let £ C Doni(t). Below we denote by Spec(r, £) the spectrum of the classical 
weak Galerkin problem: 

3 u E 'C\{0} and /i e R such that i{u,v) — iJ,{u,v) for all v E £. 

Assume that an interval (a, 6) C M is such that Tr[/^_^ dE\] = Tr[/^j^ dE\] = 
oo. Then for general (£„) E A(r), the set 



(a,6)n fl y Spec(r,£fe) 

n=l k>n 

could be a much larger set than (a, b) n Spec(T). This phenomenon is usually called 
spectral pollution. See [22] for further details on this in case T — Ms and £„ chosen 
as finite element spaces. 

The following classical convergence result will play a fundamental role below, 
see [SJ Theorem 6.11]. Assume that ||T|| < oo. Let (£„) E A. For any isolated 
eigenvalue A E Spec(r) \ conv[SpeC(,ss(T)], 

(6) ^(ker(T - A), £„) - £„ ^ dist[A, Spec(T, £„)] - 0{el). 

If T is only bounded from below, the condition £„ C Dom(t) and (£„) E A, is 
typically not sufficient to ensure approximation. By applying the spectral mapping 
theorem it can be shown that ((G]) still holds true for a A < min Spec^gg (T) whenever 
(£n) G A(t) and 5 is replaced by 5i, see for example the trick applied in the proof 
of [21 Corollary 3.6]. This type of convergence is often called superconvergence. 

3.2. The Zimmermann-Mertins method. The following method for computing 
eigenvalue enclosures originated from [19 and is closely related to the classical 
Lehmann method. Below we show that it may be described in simple terms by 
means of mapping theorems at the level of reducing spaces for the resolvent. As we 
will see subsequently, convergence estimates will follow easily from ([6]). 

This method turns out to be efhcient for computing eigenvalue enclosures for Ms 
and Mc in their original matrix formulation. In Section 2] we will discuss a further 
technique which allows improvement in accuracy for Ms and depends on re- writing 
the eigenvalue problem in terms of the Schur complement. 

According to [6] Theorem 11] the present approach is equivalent with optimal 
constant to another method formulated in j4j. The latter is closely related to the 
so-called second order relative spectrum which was applied to Ms in [M]. We 
should stress that the latter is probably best for obtaining preliminary information 
about the spectrum, [3l [25] . This a priori information includes a reliable guess on 
the interval [a, b] below. 

Let a < & be such that (a, b) H Spec(T) ^ 0. Assume that a finite-dimensional 
subspace £ C Dom(T) is such that 

. {Tu,u) , , {Tu,u) 
7) mm ,1 < and max ,, ,,„ > a. 

Note that this condition is certainly satisfied by £ = £„ for n sufficiently large, if 
(£„) E A(i). We define two inverse residuals associated to the interval (a, &): 

, ((T — a)u,u) , . ((T — b)u,u) 
(8) T+ = max — ; ' ' and ^ ' ' 



use {{T - a)w, (T - a)u) ugC ({T - b)u, (T - b)u) ' 

By virtue of ([7]), we have t+ > and t_ < 0. 
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Lemma 3.1. Suppose that [a, b] D Spec(T') — {A} and a < X < b. For any subspace 
C C Dom(r) satisfying ([7]) the inverse residuals ([8]) are such that 



(9) 



6+ — < A < a 

T 



Proof. We prove the first inequahty, the second my be proved similarly. Let C 
[T - b)C, then 

(A-6)-i =min[Spec(r-fe)-i] = min ~„^L'"'' 



< min 

vec 



{{T-b)-^v,v) 



. {{T~b)u,u) 
mm 

uec {{T-b)u,{T^b)u) 



□ 



Note that 

=min[Spcc((T-6)-\£)] and t+ = max[Spec((T - a)-\ (T - a)£)]. 

This observation turns out to be useful when studying convergence of the enclosure 
© as £ increases in dimension. Below = t"^ for ~ C. 

Lemma 3.2. Let X G Spec^ig(T) and assume that [a, 6] n Spec(T) — {A} with 
a < X < b. Let *B be an orthonormal basis ofkeilT— A), and (£„) G A(T) be such 
that distT[S, £„] 0. For all sufficiently large n £ N 



(10) 6+^<A<a+^ 



1 



— - b+— =Oiei) 



Proof. The condition ([7]) is satisfied for all sufficiently large n e N, so the left hand 
side of PH)) follows from Lemma 13.11 

Let Cn = {T - b)Cn. Since (£„) G A(T), it follows that (£„) e A. Let <8 = 
{01, . . . , Then there exist vectors u„j- e £„, such that || (r — — it„.j)|| < 
(1 + !&!)£„ for each 1 < j < ^- Set u^.j = (T — b)un j G £„, then for any normalised 
(j> S ker(T — A) we have 

k fc , , ,, 



A-5 



< 



(T^-^)E 
fc(i - 



A-6 



i=i 



A-6 



^n,j ) 



6- A 



Thus (5(ker(T— A),£„) < /c(l + |6|)£„/(6— A). Hence, applying © to the operator 
(T — 6)^^ and eigenvalue (A — b)~^ yields |r,7 — (A — — 0{e^), and therefore 
b + — X = 0{e^). Similarly we have a + ^ — A = ©(e^), and the right hand 
side o"f (HOl) follows. " □ 



By means of an example we now show that (£„) G A(t) and (£„) C Doni(T), 
is not generally sufficient to ensure a decrease in the size of the enclosure as n — > 
oo. The crucial point here is that (£„) G A(t) does not ensure that distT(*B,£„) 
decrease to 0. 
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Example 3.3. Let T = Ms be as in Example EH For j g N, let (j)3k-2 = 4>k , 
03fc-i = <Pk ^nd (f>3k = , where the right hand sides are given by (|4]). Let A3fc_2 — 
Xf, , Asfe-i = 1 and Aa^ = A^. Then T(/)j = Aj(/>j and {(j)]} is an orthonormal basis 
of [L2(0,1)]3. 

For n > 2 consider the subspaces — span{Q;„0i +£„03„, (j)2, • ■ • , 03n-i} where 
e„ = p_ and a„ = ■\/ 1 — e^j . Then £„ C Doni(T) for every n e N. We show that 
(£„) G A(t). Let u e Dom(t) and 7^ = {u,<f)j). Then 



hilt 



. E(l + ^^-)l7.f <oo. 



Let M„ = 7i(a«'/'i + Sn'/'sn) + Y.^l2 ^ ^j^j- Then 

00 

\\U-Un\\1= 7l(a« - 1)01 + (7l£n - 73n)03ri - ^ 7j'/'j 

j=3n+l 



= (l + Ar)|7i|'l«n-l|' + (l + A+)|7ie«-73n|'+ J2 (l + ^J-)l7,f 

j=3n+l 

< (1 + Ar)|7ina„ - Ip + ^^ + 5,^^'' + E (1 + A.)l7.f 0, 

therefore G A(t). 

This ensures that distt((/'i, ^ 0. On the other hand, a straightforward 
calculation shows that distrC^i, 'Cn) — > 1 + a/3 as n -> 00. Choose a = and 

b — 1^^'^ so that [a, 6] n Spec(r) = {Aj^}. For n large enough, 

^ b = a„Ai + — " ( o ) < 



la„0i + e„(?!)3„||2 " A„ 

{T(l)n-l,(t>n-l) 



- a = A+ > 0. 



Thus ([T]) is satisfied and we may obtain upper and lower bounds on the eigenvalue 
Aj^ from the left side of (ITUl) . 

Let us now prove that the length of the enclosure in ([TU|) does not decrease. 
Indeed 

aliX^r + l ^(AD^ + l^^^^^ 



as n ^ 00. 

It is easily verified that 

ST{ker{T-\),Cn)^0{e,,) ^i(kcr(r - A), £„) - 0(£„). 

As the following example shows, 5T(ker(r — A),£„) and (5t(ker(r — A),£„) can 
converge at the same rate, but the latter may be faster. Therefore, there is a 
potential loss in convergence of the method when compared with the standard 
Galerkin method in the case where the latter is applicable. This loss of convergence 
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is compensated by the fact that the enclosures found are certified and free from 
spectral pollution. 



an = 



Example 3.4. Let T and 0„ be as in Example 13.31 Let e„ — (Aj) ^ and 
y/T~£l. If we consider £„ = span{a„0i + e„03n-i, 02, ■ • ■ , 4>^n-2}, then (£„) e 
A(T), and both 5t(ker(r-l),/:„) and (5T(ker(r- 1), £„) areC'(n-2), If we consider 
Cn = span{a„0i + e„03„, (/)2, . . . , (/)3„_i}, then (£„) e A(T) once again but now 
(5i(ker(T- AJ"),£„) = 0{n-^) and ^T(kcr(r - A^), £„) = 0{n-'^). 

4. Operator matrices and eigenvalue enclosures 

The linearized MHD operator associated to the plasma layer configuration ^ 
falls into the class of top dominant block matrices. We show that enclosures for 
the eigenvalues of Ms which lie above the essential spectrum can be obtained from 
enclosures for the eigenvalues of its Schur complement. Denoted by S{^), the latter 
is a /Lt-dependant holomorphic family of semi-bounded operators. Upper bounds for 
its eigenvalues can be found from a direct application of the Galerkin method. We 
show in Section 14.21 that these upper bounds are superconvergent as the dimension 
of Cn increases, hence they turn out to be asymptotically sharper than the upper 
bounds found from the method of Section [3] applied directly to Mg. In Section [1751 
on the other hand, we show how to find lower bounds for the eigenvalues of Mg 
from corresponding lower bounds on the eigenvalues of S{^). The latter are found 
from the left side of ^ with T = S{fj.) for a particular choice of fi. 

4.1. Basic notation. The results established below apply to any block operator 
matrix Mq as in ^ which is top dominant in the following precise sense, see f26]. 

a) A and D are self-adjoint operators acting on Hilbert spaces Hi and ^2, 
respectively. 

b) Ais bounded from below, D is bounded from above, B is closed and densely 
defined on a domain of 'H2 with values in Hi- 

c) Dom(|A|2) c Dom(B*), Dom(B) C Dom(i:») and Dom(B) is a core for D. 

Without further mention, we assume that the entries of Mq are subject to these 
conditions. They are satisfied by the plane slab configuration MHD operator, how- 



ever, for TO 7^ the ansatz c) does not hold in general for the cylindrical pinch 
configuration. 



The first condition in and the semi-boundedness of A, together, imply the 



existence of constants a, /3 > such that 

(11) \\B*uf < aa[u]+ /3\\uf for all u G Dom(|^|5) = Dom(a) 

where a is the closure of the quadratic form associated to A. These two constraints 
also imply that {A — i')~^B is a bounded operator, so Dom{{A — i')~^B) = H2, 
for an arbitrary i' < minSpec(A). The self-adjoint closure of Mq, which we denote 
here by M, is explicitly given by 

(12) Doni(A//) = I ^ ^ ^ : y G Dom(i:>), x+{A-iy)-^By G Dom(^)| 

. . / X \ / Aix+{A-iy)-^By) - u{A-v)-^By 

^^^^ \y)'\ B*x + Dy 

see [II Section 4.2] and references therein. 
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Set d — maxSpec(Z3) and U — {z & C : Re z > d}. For ^ & U consider the 
following family of forms 

(14) s{^i)[x, y] = a[x, y] - ^l{x, y) - {{D - fi)-'B*x, B*y) 

with common domain Dom(s) = Dom(s(/i)) = Dom(a). Then is a holomor- 
phic family of type (a), see JJ:, Proposition 2.2]. Associated to these forms is a 
holomorphic family of type (B) sectorial operators 5'(/i): 



(15) Dom(5(Ai)) = {a; e Dom(a) : x - (A - v)-^B{D - ^jlY^B*x E Dom(A)} 



(16) S{^l)x = {A- v){x -{A- v)-^B(D - n)-^B*x) + {v - ^i)x] 

see [71 Proposition 4.4]. Here, as above, v < min Spec(A) is fixed, but can be chosen 
arbitrarilly. We note that for any x e Dom(s) and fi G J7 n R, 

(17) ^£|M^_||^||2_||(^_^)-l5*^||2, 

The families s(-) and 5(-) are called the Schur form and the Schur complement 
associated to Mq, respectively. 

The form s(/i) is symmetric and semi-bounded whenever fi E M.nU. The corre- 
sponding operator S{fi) is therefore self-adjoint and bounded from below. We set 
the spectra of the Schur complement as 

Spec(5) = {fieU :0e Spcc(5(Ai))}, 
SpeCdis(5) ^{tJ.eU :0e Spec^i,{S{ii))} and 
Spec,,,(5) = {^IeU■.Oe Spec,,,{S{fi))}. 

4.2. Upper bounds via Schur complement. We denote 

Ac = inf{SpeCcss(M) n (d,oo)} 

and Ai < A2 < • • • the repeated eigenvalues of M which lie in the interval {d, Ac). 

Lemma 4.1. The spectra of S and M coincide on (d, Ae) and SpeCgss(5')n (d, Ac) — 
0. Moreover dimKer(S'(Aj)) = dimKcr(A/ — Aj). 

Proof. For the first and third assertions, see Proposition 4.4]. For the second 
assertion we proceed by contradiction. 

Suppose there exists A G {d, Ac) such that £ SpeCcss(S'(A)). Since S{X) — S{X)* , 
there is a singular Weyl sequence a;„ G Dom(S'(A)) such that ||a;„|| = 1, Xn 
and ||S'(A)a;„|| — ^ 0. Let 

{-{D-Xi)-^B*Xr. 
Then y„ e Dom(Af) and ||j/„|| > 1. A direct calculation shows that 

(M-A)2/„= (^^(^q)^") ^0. 

We prove that 2/„ has a subsequence yn{k) ~^ 0, which in turn is a contradiction 
because A ^ SpeCcss(M). Let V = {D — A)Dom(i?). By virtue of the second and 
third ansatz in 2? is a dense subspace of T-L2- Moreover, 

(18) {{D - \)-^B*Xn,y) -^Q for all ye P. 
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According to ([TT1) . 

\\B*Xn\\'^ < aa[xn] + 13 < a(o[x„] ~ X) + aX + l3 

< a5{X)[xn] + aX + /3 — a{S{X)xn, Xn) + aA + /3. 

As the right hand side of this identity is uniformly bounded for all n, there exists a 
subsequence and z £ 'H2 such that {D — X)^^ B* Xn(^k) z. Since (fT8|l implies 
that z = 0, the subsequence y„(fc) is as needed. □ 

For pi > d, we denote the spectral subspace of S{^) corresponding to an interval 
J by 



£j{S{fi)) ^ Range J^dEx 



Here we abuse the notation and write dEx for the spectral measure associated to 
the self-adjoint operator S{^) also. Let the dimension of £(_oo,o) ("S'Im)) be 



k{^) — tr dE\. 

J(-oo.O) 



'(-00,0) 

Throughout this section we assume that < 00 for some fi > d. By [71 Theorem 
4.5], this assumption and aj|[cj imply the existence of 7 > d such that (d, 7] n 
Spec(Af ) = 0. We write k := ^(7) < 00 and note that k is independent of the 
particular choice of 7 € (d, min{Spec(M) n (d, 00)}). 

Let li{m) — min{j N : Xj — Am} and him) ~ max{j Cz N : Xj — A„i}- Then 
k(-) is constant on intervals contained in (d, Ac) \ Spec(Af ) and 

(19) dimf(_oo_o)('5'(Am)) = k + li{m) - 1, 

(20) dim£(_oo,o]('S'(Am)) = K + l2{m); 
see [II Section 2] for further details. 

Example 4.2. In the case of the plasma layer configuration, S{fj,) is a family of 
Sturm-Liouville operators. It is readily seen from the results of [TTl §7.5] that 
Ac — 00, K < 00 and Spec(M) n (d, Ae) consists of a sequence of simple eigenvalues 
which accumulate at +00. 

We now describe the theoretical framework and basic procedure for approximat- 
ing a fixed eigenvalue Am • Denote by Ei (p.) < • • • < i?K+m the first k + m 
eigenvalues of S{fi) repeated according to their multiplicity. Let 

C — spanjui, . . . , Un} C Dom(s) where {ui, uj) = 5ij, 

be an n-dimensional subspace where n > k -\- m. Consider the family of matrices 
ScilA ^ C"^" whose entries are given by 

Sc{fJ')ij = s(p)[uj,Ui\ i,j = l,...,n. 

Denote by Ei{£,,fi) < ■■■ < E^+mi^, fJ-) the first k + m eigenvalues of ScifJ-) 
repeated according to their multiplicity. 

Lemma 4.3. Let 11 e (d, 00) be such that iJ^+ml'C, fi) < 0, then Am < M- 

Proof. We suppose that Am > The Rayleigh-Ritz variational principle ensures 
that Eii^mdj-) < -EK+mC-C, /i) < 0, and therefore 

(21) dim£(^_ao,o]{S{tJ')) > K + m. 
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If /i = Xj for some some 1 < j < h{m), then from ([20|) we have duTi£(_oo,o]('S'(A*)) = 
K+hU) < K+li{m), which contradicts (|21|) . Suppose now that /i ^ Spec(M). Then 
from ([T7| we have s(/i) > s(Am) from which it follows that dim£(_oo o)(*S'(/i)) < 
dim£'(_oo,o)('S'(Am))- From (fT9|) we then deduce that dmi£(^_oo,o){S{^)) < k + m, 



An upper bound for Am may be obtained by applying the Galerkin method to 
the Schur complement, then finding a G (d, oo) such that S{fi) has at least n + m 
non-positive eigenvalues via a root finding algorithm. We now turn our attention 
to the convergence properties of this approach. For this we employ ([6]) assuming 
T = 5(A,„) and denote £ ^ ker(r). 

Theorem 4.4. Let (Cn) G A(t) be such that 5i{£,Cn) = En ^ ^ as n oo. Let 

/i^ G M &e such that Ei^^mi^n, f-n) = 0- Then /i+ > A„j and /i+ — A„i = ©(e^J. 

Proof. From ([T9|) it follows that S{Xm) has k + li{m) — 1 negative eigenvalues 
counting multiplicity. Therefore, the density condition (£„) S A(t) implies that 
for all sufficiently large n G N there are precisely k + /i(m) — 1 elements from 
Spec(T, £„) which are negative. Since di{£,Cn) 0, there are precisely 

him) — hifn) + 1 (= dimf) elements from Spec(T, £„) which are non-negative and 
of the order ©(e^). The result now follows from ([6]) and ([TTl) . □ 

4.3. Lower bounds via Schur complement. In the previous section we found 
upper bounds for an eigenvalue A„i G (d, Ae). We now turn our attention to find- 
ing complementary lower bounds for this eigenvalue via the method described in 
Section [321 

Lemma 4.5. For ^ G {d, Ac) let < e < ^ - d. If [-e, 0] n Spec(S'(/i)) 7^ 0, then 
[fi - n Spec(M) ^ 0. 

Proof. Let 6 = max{Spec(S'(/i))n[— e, 0]} and assume that [fi—e, /i]nSpec(M) = 0. 
According to [H Lemma 2.6], k(^) = — e). By virtue of (|T7|) and the Rayleigh- 
Ritz variational principle, 

> mf sup ,, >e+ mf sup , „ £ + (5, 

where the right hand side is non-negative. The result follows from the contradiction. 



By applying Lemma 14.31 we can find /i G K such that X„i < fJ,. If S{pi) has 
K + him) non-positive eigenvalues and 6 > is such that (0, h] n Spec(5'(/i)) = 0, 
then we employ the Zimmermann-Mertins method with T = to obtain a lower 
bound on the first non-positive eigenvalue. Combined with Lemma 14.51 this yields 
a lower bound for A^. We now find the rate of convergence of this lower bound. 

Lemma 4.6. Lei > 6e smc/i </iat (0, 6] nSpec(5(A„j)) = 0. Let he a sequence 
of real numbers such that < /i„ — A™ = as n —7> 00 . For all n sufficiently 

large (0,6] n Spec(S'(/i„)) = 0. Moreover, 



which contradicts ([21]) . 



□ 



□ 



(22) 



|1(5(A,„) - 6)-i - iSi^in) - h)-^\\ = 0(£„). 
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Proof. Wc first show that b ^ Spec(S'(/i„)) for aU sufBciently large rt g N, and that 
([22l) holds true. Let x E Dom(s) and a„ = /!„ — Am, then 

s{\„,)[x] = a[x] - X,n\\x\\^ - {{D - X„,)-^B*x,B*x) 

= s(/i„)H + an\\xf + {[{D - Ai„)-i -{D- Xrn)-^]B*x, B*x) 

= s{iin)[x] + a„||x||2 + - Mn)"^^ - A„0"iB*x, B*x). 

Set s — s{fin) — s(Am). Note that a[x] < s{X)[x] + X\\x\\^. By virtue of (fTTI) . we 
have 

|s[x]| = a„||a;||2 + a„((i? - ^„)-i(i^ - Xmr'B*x,B*x) 

< ajxf + ar,\\{D - fi^yHD - Xrny'\\\\B*xf 

< an\\xf + aJ{D - tiny\D - X„^y^\\{aa[x] + /3||xf) 

< a„(a„||x|p + bnS{Xm)[x]). 



where 



1 + (3\\{D - f^ny\D - A™)-i|| + aA„ ^ 1 + /3(A™ - d)-^ + aA„ 



bn ^ a\\{D - ^in) ^{D-Xrn) ^||^a(A„-d)" 



as n — >■ cx). 



Set ci = niax{||(S'(Am) - by^\\,\\S{X,n){S{X^) - by^\\}, and let n g N be sufB- 
ciently large to ensure that a„(a„ + &„) < c^^. By virtue of [T21 Theorem VI-3.9], 
we obtain b ^ Spec(S'(/i„)) and 

(23) \\{S{X„,)-by'-{S{^l^)-by'\\ < 4a„(a„ + Mc? 



(1 - a„(a„ + 6„)ci)2 

which immediately implies (1^^ . 

It remains to show that (0, b] D Spec(5(/i„)) = for all sufficiently large n E N. 
By virtue of [121 Theorem VII-4.2], there exists a constant < min Spec(5(/iri)) 
for all sufficiently large n G N. Let F be a circle with center (& + i/)/2 and radius 
{b — i')/2, and set 



C2 — max 



{|i(5(A,„) ~ ||5(A„0(5(A,„) - 



Then an{an+bn) < C2 ^ for all sufficiently large n g N. Applying the same argument 
as above, we obtain 

max^ ||(5(A„J - z)- - - z)-^\\ ^ < + ^"^'^^ 



;{||(5(A„0-z)-1-(5(a*„)-^)-'I|} 



The right hand side of converges to zero as n — > 00. Thus, the spectral subspaces of 
S(Xrii) and S(fin) corresponding to the eigenvalues below b have the same dimension 
for all sufficently large rt g N. The desired conclusion follows from ([TT)) and the 
Rayleigh-Ritz variational principle. □ 

According to Theorem l44l if (£„) g A(s(A„i)) and Ss(x^)i£,Cn) = 0{en), then 
we obtain a sequence of upper bounds 

(24) fi+\Xm satisfying ^+ - A„ 0(4)- 



If we now pick b satisfying the hypothesis of Lemma 14.61 by virtue of (|9]) we find 
a lower bound for the smallest in modulus non-positive eigenvalue of T = S'(At^). 
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That is 

6 + — < max {Spec(S'(^+)) n (-00, 0]} 

Tn 



Lemma 14.51 ensures corresponding lower bounds 

(25) = ii+ + b+l/T- <\n. 

In the theorem below we find bounds on the speed of convergence /i~ Xm- 

Before proceeding further, we note that S{ii) and B{D — Xm)~^(D — fj,)~^B* 
for the plane slab configuration ([2]) are sectorial Sturm-Liouville operators for all 
fi E U . Both families of operators are closed in the domain 

(26) V = H\{Q, 1); dx) n HI{{Q, 1); dx), 

which coincides with P3|) and is independent of Moreover, they are both holo- 
morphic families of type (A), see [121 Example Vn-2.12]. 

Theorem 4.7. Suppose that the entries of Mq satisfy^^ -[c)| Assume that S{^) is 
a holomorphic family of type (A) with Dom(S'(/z)) = V independent of fi. Assume 
additionally that B{D — Am)^^(-D — fi)-^B* is closed on D. Let (Cn) G A{S{Xm)) 
with 5s(\^){£,Cn) = Sn- If l^n Constructed as in ((25l) . then fi^ < Xm and 
Xm - IJ-n = C'(e^) as 00. 

Proof. Let /i G (d, Ac), x £ V and ly < minSpec(^). According to ^TE\i we have 
S{^i)x - SiXm)x = (A - iy){x - {A-iy)-^B{D - n)-^B*x) + {v - ti)x 

- (A - iy){x - {A-u)-^B{D - Xmr^B*x) -{v~ Xm)x 



= (A - v){fA - v)-^B\(p - Xm)-^ -{D- ^i)-^]B*x) 

+ (A„i - ^i)x 

= (Xm - ^i){A - u){{A-v)-^B{D - Xm)'\D - y)-^B*x) 
+ (A„i - /i)a; 

= (A,„ - n)B{D - Xm)-\D - ^i)-^B*x + {Xm - ii)x. 
We consider the closed operator B{D — Xm)-^{D — fj,)-^B*. Since 

B{D-Xm)-HD-f,)-^B*x = ^^l^^^^^^^-x for ^ ^ A„„ 

Am — M 

and S{fi) is holomorphic of type (A), then \\B{D — Xm)-^ (D — tJ-)-^ Bx\\ is uniformly 
bounded in a neighbourhood of /i. Moreover, for any y e Dom(_B*) the function 

{B{D - Xm)-\D - ^i)-'B*x, y) = {{D - ^,)-'B*x, {D - Xmy^B*y) 

is analytic for jj. E U. It follows that B{D — Xm)-^{D — nj^^B* is a holomorphic 
family of type (A) . 

Let J C (d, Ac) be any compact interval containing a neighbourhood of Am- By 
virtue of [T2| Section VII. 2.1], there always exists a constant C3 > 1 such that 

(07^ \\x\\ + \\B{D-Xm)-HD-l^r'B*x\\ 

(27) — — -— — - — ^ < C3 tor all u E J. 

^ ' \\x\\ + \\B{D-Xm)-^B*x\\ - ^ 

Since the operators B{D — Xm)-^B* and S'(Am) have the same domain V, there 
exist constants a, /3 > such that 

(28) \\B{D -Xmy^B*x\\<a\\S{Xm)x\\+P\\A\ for all x eV. 
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As ^+ — A = C(£^i), for some N £ N large enough e J whenever n > N. 
Combining ^ with ([27]), gives 

\\B{D - X,n)-\D - fi+)-'B*x\\ < (C3 - 1)11x11 + C3\\BiD - Xra)-^B*x\\ 

< (c3 - l)||a;|| + ac3\\S{Xm)x\\ + f3c3\\x\\ 

<C4(||5(A,„)a;|| +||.t||), 

where C4 > is independent of n > A^. Thus 

(29) \\S{^l+)x-S{XM<{^^n~^r,^){cA + lmS{XM + M) forall n>N. 

Let Cn — (S'(/x+) — b)Cn- We show that (£„) G A. Let v £ H. There exists 
u & V such that {S{Xm) — h)u = v. Since (£„) g A(5(A„i)) we have a sequence 
Un £ £n satisfying {S{Xm) — b)un — > w. As b ^ Spec(5(Am)), the sequences 
and ||S'(Am)M„|| are uniformly bounded. Hence, it follows from (|29p that 
\\S{p+)un - 5'(A„)u„|| 0. Since 

V - (S'(^+) - b)Un + {S{^J.'^) - S{X,n))Un = V- (S'(A,„) - fe)w„ -> 0, 

clearly also (5(^+) — 6)u„ w. Thus £„ e A. 

We now show that S{£,Cn) — C(£n)- Let ipi, . . . , ipk be an orthonormal basis 
for £. There exist vectors Unj G Cm such that \\S{Xm){4>j — Unj)\\ < and 
\\4>j ~ Wnjll < e„ for each 1 < j < fc. We set Unj = {S{^n) ^ b)unj e £„. Using 



291) we have for any normalised <j) € £ 

k 



k 



< 



k 



E 



^(5(M+)-&)u„, 



+IIE' 



3 = 1 





1 

^E(ii'^('^«)^"j~ 


5'(A™)u„j|| + \\S{Xm)Un,j\\) 


IkSfi 


1 

+ ^E(^" - A,„)(C4 
i=i 


+ 1)(|15(A„)m„,,|| + ||w„,,|1) 




1 


+ !)(£„ + ||ltnj||)- 



Therefore (5(£, £„) = C'(e„). 

We complete the proof of the theorem as follows. By applying ^ to the operator 
T = {S{Xm) ~ b)^^ and eigenvalue (— 6)~^ — min{Spec((S'(A„i) — 6)^^)}, we obtain 

(30) min{Spec((5(A™) - &)-\£„)} + - 0(4). 
Using Lemma HTEl and < /i+ — A™ = ©(e^), we have 

(31) \\{S{X^) ~ b)-' - (5(M+) - &)-^|| = 0{el). 

From (1501) . (PT|) and the Rayleigh-Ritz variational principle, it becomes clear that 

(32) T- = min{Spec((5(^+) - 6)"\£„)} satisfies +6-^ = 0{el). 
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Figure 1. Log-log graph. Vertical axis: 



X 100. Horizon- 



tal axis: eigenvalue index m. We depict the relative size of the 
enclosure in the calculation of the first 50 eigenvalues of Mg of Ex- 
aniple[121 The subspace £„ is chosen to be: C^f and C{h, 1, r) for 
Hermite elements of order r = 3, 4 and 5 on an uniform mesh with 
h chosen so the dimension of the spaces is approximately 10 x 67. 
We have chosen n = 10^^"^, Ts = 10^^^, Tb ~ 0(10^^) in the case 
of the sine basis, and n = 10"^", Ts = 10"^ and Tb = 0(10"^) in 
the case of the finite element method. 



Moreover, b + 1/t~ is precisely the lower bound on the smallest in modulus non- 
positive eigenvalue of S{fi^) which is obtained from the Zimmermann-Mertins 
method. Then /i~ = + b+l/T~ < follows from Lemma and X„i — ^i" = 
0{el) follows from 1^. □ 



5. Numerical examples: plane slab configuration 

An optimal strategy in terms of convergence for calculating enclosures for the 
configuration ([2]) can be established from the approach in Section |4l We now 
illustrate the practical applicability of this strategy by performing various numerical 
experiments on benchmark models. Our equilibrium quantities will be chosen from 
examples 12.11 and 12.21 

For a fixed ji e [d, oo), the eigenvalues of S{^) are simple in both examples. The 
corresponding eigenvectors are in C°°(0, 1) and they satisfy Dirichlet boundary 
conditions at the endpoints of the interval. We have 



Spec(Ms) n (d, oo) = {Ai < A2 < . . .} 
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Max element size 



Figure 2. Log-log graph. Vertical axis: \^i23~ h\- Horizontal 
axis: maximum element size h. The subspace is chosen to 
be: 



1,3) for decreasing values of h. For these calculations 



10 



-10 



10 



-12 



Tb = 0(10-3). 



where each eigenvalue is simple and Xj —>■ oo. Below we distinguish the model 
used by denoting these eigenvalues by and ASP' respectively. For all u S 

Dom(s) = H^iiO,l);dx), 



(33) 



> vr^llulp + ((7/4 - x)u, u) ~ fi\\u\\ 

> {n' - ^l)\\u\\^ 



Thus S{^) is positive definite for /i £ (d, tt^). Hence upper bounds /i^ for Am can 
be found from Theorem 14.41 with k = 0. 

In practice we find /x+ as follows. For a fixed £„ we compute a few eigenvalues 
of S{fi) for n in an uniform partition with p points of a suitable interval (a, b) 
containing only Am- We then approximate /i^ via one iteration of Newton's method. 
Below, the integrations involved in the assembling of the matrix problems are set 
to a tolerance of the order ©(tj), the eigenvalue solver is set to a tolerance of order 
0{ts) and Tb = (6 — a) /p. These are different for the different experiments. By 
virtue of pT|) . the root finding step is accurate to 0{t^). 

To find complementary lower bounds from Theorem 14.71 we require b > such 
that (0, b] n Spec(5'(/x+)) = for all sufficiently large n e N. From dSS]) it follows 
that s(^+)[m] > a[u\ — fi'^\\u\\^ for all u £ Dom(s) (in both Example 12.11 and 
Example 12. 2p . By the minmax principle, the (to -|- l)-th eigenvalue of S'(^+) -I- /i+ 
lies above the (to-I- l)-th eigenvalue of A. In fact A™ < {m + l)'^TT'^ < Xm+i and we 
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may choose b E (0, (m + l)^7r^ — Integration and eigenvalue solver tolerances 
are set as for the upper bounds. 

We consider two canonical basis to generate £„ C V, see (|26p . A first natural 
choice is the sine basis, 

= spanjui, . . . , w„} where u„ = V2sin{mTx). 

Standard arguments show that £„ e A{S{Xm)) and 

where r can be chosen arbitrarilly large. Applying Theorem 14.41 and Theorem 14.71 
we obtain 

\ Am, A^n/^Am, and /i+ - /^^ = 0(n"''). 

This means that the enclosures should converge to zero super-polynomially fast for 
the family of subspaces See Tableland Figure [TJ All calculations involving 

this basis were coded in Matlab. 



n 








A4 


^5 


5 


1/;.OOU38099 


41.9106P™ 


91.2474™fi 


160.3305?!?^ 


249.155j?og8 


10 


12.3504^11 


41.910611?^ 


91.2474™3i 


160.33058?!! 


249.15507911 


20 


12.3504lili 


41.91061214 


91.24747026 


160.330582[j| 


249.155079|i 


40 


12.3504|i^ 


41.91065^J? 


91.2474™26 


160.330582|} 


249.155079^^ 



Table 1. Approximation of the first five eigenvalues of Afg for 
Example 12.21 and test spaces chosen as C^^ . For these calculations 

Ti = 10-14, T, = 10-12, Tb = 0(10-5). 



Another natural basis is obtained by applying the finite element method. Let 
S be an equidistant partition of [0, 1] into n sub-intervals /; = of length 

h = l/n = xi — xi^i. Consider the subspaces 
(34) 

Cih, k, r) = Vk{k, r, E)^{ve C'=(0, 1) Pr{Ii), 1 < I < n, i;(0) = = i;(l)}. 

The C{h, k, r) are the finite element spaces generated by C'^-conforming elements 
of order r subject to Dirichlet boundary conditions at and 1. Then 

\\V - VhWwP.^O,!) < C||w|| Hg+i(o,l)^''^^"^ 

for Vh e jO.{h, k, r) the finite element interpolant of u G C*"' n H~'^^{0, 1). For fixed 
k,r > 1, let /i^j ^ and ^ be the upper and lower bounds for Am given by Theorem 
14.41 and Theorem 14.71 respectively. Then 

(35) - Am - 0(/i''') and f,- ,^-X^=Oih^^^-'^). 

All calculations involving this basis were coded in Comsol. 

Figure [5] shows that the orders of convergence found in are optimal in the 
case r = 3 for Hermite elements and to = 1,2,3. In order to compare the quality 
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Eigenvalue index 



Figure 3. Log-log graph. Vertical axis: ^" x 100. Horizon- 
tal axis: eigenvalue index m. We depict the relative size of the 
enclosure in the calculation of the first 50 eigenvalues of Mg of 
Example 12.21 The subspace £„ is chosen to be: £.{h, 1,5) for de- 
creasing values of h. For these calculations Tj = 10^^", Tg = 10~^, 
^ = 0(10-3). 



of the upper and lower bounds, we have chosen Example 12.11 and calculated the 
value of with the exact formula in machine precision. Observe that the upper 
bounds are all roughly 4 orders of magnitude more accurate than the lower bounds. 
This is certainly expected from the fact that the calculation of the upper bound 
involves the solution of a second order problem, whereas that of the lower bound 
involves the eigenproblem (|8]) with T — S{fi^) which is of fourth order. Here we 
have purposely chosen large values of h, so the calculation of the bounds for A2 and 
A3 is not particularly accurate. 

The aim of the experiment performed in Figure [T] is to compare accuracies in 
the computation of the bounds by picking of roughly the same dimension, but 
generated by different bases. For this we have fixed £„ of a given dimension and 
compute the size of the enclosure relative to the size of the lower bound 

fi~. We consider Example 12.21 We have chosen dim£„ = 67 for the sine basis 
and dim£„ 670 for the finite element bases (remember that the sine basis is 
exponentially accurate). 

The accuracy deteriorates (even in relative terms) as the eigenvalue counting 
number m increases. For the same dimension of £„, accuracy increases as the 
order of the polynomial r increases. In this figure, the enclosures found for Am for 
TO < 5 (r = 3), m < 10 (r = 4) and m < 20 (r = 5) should not be trusted and it 
is just included for illustration purposes. This locking effect is consistent with the 
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fact that the calculation of the enclosures can never be more accurate than a factor 
of max{ri,Ts,T^}. 

We can examine this phenomenon in more detail from Figure [3] and the blue line 
in Figure [2l As the dimension of the test subspace decreases, for each individual 
eigenvalue, the residual starts decreasing and eventually hits the accuracy threshold. 
From Figure [2] it should be noted that the lower bound hits the threshold earlier 
than the upper bound, however this threshold for the lower bound is three to four 
orders of magnitude larger that that of the lower bound. 



6. Numerical examples: cylindrical pinch configuration 

The approach considered in Section |4] cannot be implemented on the cylindri- 
cal pinch configuration for m as the block operator matrix does not satisfy 



condition c ) We now report on a set of numerical experiments performed on the 
benchmark model in Example 12. 3[ by directly applying the method described in 
Section El to T = M. 

In this case we have chosen £„ = C{h, 1, r) x C{h, 1, r) where C{h, 1, r) is defined 
by (I34p and is generated by Hermite elements. The Dirichlet boundary condition 
inposed at both ends of the interval [0, 1] ensures that £„ G Dom(Af). In Tabled 
we show computation of the first three eigenvalues above SpeCpj,g(M). Similar 
calculations can be found in }101 Table 1]. Note that in the latter, for iV = 32 the 
approximated eigenvalue appears to be below Ai whereas for = 64 it appears to 
be above Ai. This phenomenon is not present in the method described in Section [3] 
as it always provide a certified enclosure for the eigenvalue. 



3 


exact Am 


(a,fo) 


enclosure 


d.O.f. 


1 
2 
3 


4.38995771667 
29.4242820473 
73.8686971063 


(3,20) 
(20,60) 
(60, 100) 


A q9039t)2 
^■■^895445 

29.42i56 

70 SR9378 

(0.0DgQ3Q 


5004 
5720 
8004 



Table 2. Enclosures for the first three eigenvalues in SpeC(jjs(Af) 
above the essential spectrum for Example ()2.3|) by direct applica- 
tion of the method of Section [3l For these calculations we have 
chosen r = 3 and — — 10~^. 



The eigenfunctions of M associated to Am possess a singularity at the origin, so 
neither the upper nor the lower bounds obey an estimate analogous to that of ([35]) . 
On the left of Figure 2] we show a log- log plot of the size of the enclosure against 
maximum element size for r = 3 and r = 5. The graph clearly indicates that the 
order of decrease of the enclosure does not seem to decrease with the order of the 
polynomial. On the right of Figure U we show the absolute residuals for lower and 
upper bounds separately. Both graphs indicate that 

=0(/ii) as/i^O 

equally for r = 3 and r = 5. 
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10"^ 10^ 10'^ 

Max element size Max element size 



Figure 4. Log-log graphs. Vertical axis: /i^ ^ — fi^ ^ (left) and 
|Am"/^^/j| (right). Horizontal axis: Maximum element size /i. For 
the subspaces £„ we choose r — 3,5 (left) and r = 3 (right). Note 
that the order of decrease of all the residuals is roughly 0{h^) for 
both polynomial orders (left) and both bounds (right). For these 
calculations r, = Tg = 10^^". 
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